Direct Imaging in the Habitable Zone and the Problem of Orbital Motion 

Jared R. Males, Andrew J. Skemer, Laird M. Close 
Steward Observatory, University of Arizona, Tucson, AZ 85721 
^^ I jrmales@as.arizona.edu 

ABSTRACT 



o 

^H, High contrast imaging searches for exoplanets have been conducted on 2.4-10 m telescopes, 

typically at H band (1.6/ini) and used exposure times of ~ Ihr to search for planets with semi- 
C\J ■ major axes of > 10 AU. We are beginning to plan for surveys using extreme-AO systems on 

^^ - the next generation of 30-mctcr class telescopes, where we hope to begin probing the habitable 

zones (HZs) of nearby stars. Here we highlight a heretofore ignorable problem in direct imaging: 
planets orbit their stars. Under the parameters of current surveys, orbital motion is negligible 
over the duration of a typical observation. However, this motion is not negligible when using 
^S^ ' large diameter telescopes to observe at relatively close stellar distances (1-lOpc), over the long 

Mh, exposure times (10-20 hrs) necessary for direct detection of older planets in the HZ. We show that 

Q I this motion will limit our achievable signal-to-noise ratio and degrade observational completeness. 

+j ' Even on current 8m class telescopes, orbital motion will need to be accounted for in an attempt 

^ . to detect HZ planets around the nearest sun-like stars a Cen A&B, a binary system now known 

to harbor at least one planet. Here we derive some basic tools for analyzing this problem, and 
ultimately show that the prospects are good for de-orbiting a series of shorter exposures to correct 
^ . for orbital motion. 
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VO ' 1. Introduction 

^' 

f^ . Orbital motion has been used in one fashion or another to detect planets around stars other than our 

f^ I Sun in large numbers. The radial velocity (RV) technique monitors the Doppler shift of a stellar spectrum 

as the star itself orbits the planet-star center of mass, thus allowing us to infer the presence of a planet. 
Similarly, the astrometry technique monitors the motion of the star on the sky and likewise infers the presence 

S^ I of a planet. The transit technique monitors the reduction in brightness of the star as the orbiting planet 

H ■ temporarily crosses the line of sight between the telescope and the star. 

Unlike these indirect te chniques, direct imaging dete cts light from the planet itself and spatially resolves 



rig 



it from the light of the star (jTraub fc Oppenheimeiil2011l ) . The extreme difference in brightness between star 



and planet at small projected separations has generally limited direct imaging efforts to wide separations 
where orbital motion is ignorable. The next generation of large telescopes will move us into a new regime of 
direct imaging, moving closer to the star. We will even be able to begin probing the liquid water habitable 
zone (HZ). Here we point out that at these tight separations orbital motion will no longer be negligible in 
direct imaging. As we will show the motion of planets in the HZ (and closer), during the required integration 
times, will be large enough to limit our sensitivity unless we take action to correct it. 

In Section [2] we present our motivation for this study and briefly review some of the related prior work. 
In Section [3] we develop the basic tools needed to analyze this problem, including the expected speed of 
orbital motion in the focal plane and the effect it has on signal-to- noise ratio {SNR) . In Section |4] we 



analyze the impact orbital motion will have on a search of a Cen A by the Giant Magellan Telescope (GMT) 
working at 10/ini, and propose a method to mitigate this impact by de-orbiting a sequence of observations. 
Then in Section [5] we treat the more favorable case of a cued search, where we have prior information from 
an RV detection. To do so we analyze the case of the potentially habitable planet Gl 581d being observed by 
the planned European Extremely Large Telescope (E-ELT). Finally, in Section [6l we present our conclusions 
and prospects for future work. 



2. Motivation and Related Work 

Moving the hunt for exoplanets into the HZ of nearby stars marks a departure from prior efforts. Here 
we briefly discuss the definition of the HZ, review direct imaging results to date, discuss the differences 
between them and and future efforts, and finally review some closely related prior work. 



2.1. Nearby Habitable Zones 

The HZ is generally agreed to be the region around a star where a planet can have liquid water on 
its surface. This is far from simply related to the blackbody equil ibrium temperature , as it depends on 



atmo spheric composition and the action of the greenhouse effect (JKasting et al.l Il993l : iKopparapu et al 



20131 ). among other factors. For our purposes it is enough to assume that the HZ is generally located at 



about one AU from a star, scaled by the star's luminosity 



aaz ~ VL * /Lq AU. (1) 



Traubl (|2012l ) provided three widths for the HZ based on various considerations, and then used the first 
136 days of data from the Kepler mission to estimate that the fraction of sun-like stars (spectral types FGK) 
with an earth-like planet in the HZ is r]^ sa 0.34. More generally, this analysis indicates that rjpianet ~ 1-2, 
implying that every sun-like star is likely to have a planet in its HZ, and some will have more than one. 
While this exciting result is based on a very large extrapolation from the earliest Kepler results, it is currently 
one of our best estimates of planet frequency in the HZ. 



This topic was recently brought to the fore with the announcement of a Cen Bb by iDumusaue et al. 



(|2012l ). Discovered using the RV technique, a Cen Bb is an to sin z = 1.13M0 planet orbiting a Kl star at 
0.04 AU. While certainly not in the HZ, this discovery has exciting implications for the presence of planets 
in the HZ of the nearest two sun-like stars. 

The above arguments hint that planets will be common in the HZ of sun-like stars. We are about to enter 
a new era of exoplanet direct imaging. With the next generation of giant telescopes and high-performance 
spaced-based coronagraphs we will be searching for planets in this scientifically important region around 
nearby stars. 



2.2. A DifTerent Regime 

The typical search for exoplanets with direct imaging has used 2.4m (Hubble Space Telescope, HST) to 
10m (Keck) telescopes. These surveys have mostly concentrated on young giant planets, which are expected 



to be self-luminous as they dissipate heat from their formation. This allows them to be detected at wider 
separations from their host stars, where reflected starlight would be too faint. This has also caused planet 
searches to typically work a t H band (~ 1. 6 /im), with exposure times of ^ 1 hr. Examples conforming 
to these stereotypes in clude iLowrance et al.l ( 2005[) using HST/NICMOS; the Gemini Deep Planet Search 
(ILafreniere et al. 2007 ): the Simultaneous Differential Imaging survey using the Very Large Telescope and 



MMT (|Biller et al.ll2007l ): the Lyot Project a t the Advanced El ectro-Optical System telescope (JLeconte et al 



20101 ) ; the Interna tional Deep Plan et Survey ( Vigan et al.l2012 ) ; and the Near Infrared Coronagraphic Imag' 
at Gemini South (JLiu et al.ll2010[ ) 



:er 



These searches have h ad some success. Examples include the 4 planets orbiting the A5V star HR 8799 
(JMarois et al.l l2008bl I2010I ). with projected separations of 68, 38, 24, and ~ 15 AU. These correspond to 
orbital periods of ^ 4 60, ~ 190, ~ 100, and ~ 50 years, respectively. The A5V star /3 Pic also has a planet 
( Lagrange et al.ll2010l ) orbiting at ~ 8.5 AU with a period of ~ 20 years ( Chauvin et al. 20121 ). Another A 



20081). At these wide 



star, Fomalhaut, has a candidate planet on an 872 year (115 AU) orbit (|Kalas et al 
separations it takes months, or even years, to notice orbital motion. 

In the much closer HZ, however, orbital periods will be on the order of one year. We show in some 
detail that this is fast enough to yield projected motions of significant fractions of the point spread function 
(PSF) full width at half maximum (FWHM) over the course of an integration. The resulting smeared out 
image of the planet will have a lower SNR, making our observations less sensitive. 



2.3. Long Integration Times 

In addition to HZ planets having higher orbital speeds than the current generation of imaged exoplanets, 
integration times required to detect them will be much longer. Direct imaging surveys to date have mostly 
worked in the infrared while attempting to detect young planets still cooling after formation. The coming 
campaigns to image planets in the HZ of nearby stars will focus on older planets, which will be less luminous 
in the near infrared. In the HZ, starlight reflected from the planet will be more important. The result is 
integration times required to detect such planets will be tens of hours, rather than the '^ 1 hour characteristic 
of current campaigns. 



Co nsider the Exoplanet Imaging Camera and Spectrograph (EPICS), an instrument proposed for the 



E-ELT. 



Kasper et al.l ()2010f ) predicted that EPICS will be able to i mage the RV detected planet Gl 581d , 
which has a semi-major axis of 0.22 AU with a period of ^ 67 days (IForveille et al.ll201ll:lyogt et al.ll2012l ) 



This orbit places it on the outer edge of the HZ of its M2.5V star ( von Braun et al.ll2011). EPICS will b e 
able to detect Gl 581d, at a planet/star contrast of 2.5 x 10^^, in 20 hrs with SNR = 5 ( Kasper et al.ll2010l) . 
Since this is a ground based instrument, a 20 hour integration will be broken up over at least 2 nights. 
Plausible observing scenarios could extend this to several nights, taking into account such things as the need 
for sky rotation. As we will show, the planet will move several FWHM on the EPICS detector during a 
multi-day observation. 



More generally. ICavarroc et al.l ( 20061 ) showed that when realistic non-common path wavefront errors are 
taken into account, the integration times required to achieve the 10~^ to 10^^" contrast necessary to detect 
an earth-like planet around a sun-like star approach 100 hours on the ground, even on a 100m telescope 
with extreme- AO and a perfect coronagraph. One of several concerns about the feasibility of a 100 hour 
observation from the ground is that such a long observation will be broken up over many nights. 



With net exposure times of 20 to 100 hrs, and total elapsed times for ground based observations of 
several to tens of days, HZ planets will move significantly over the course of a detection attempt. The focus 
of this investigation is the impact of the orbital motion of a potentially detectable planet on sensitivity. 



2.4. Related Work 

Though it has not yet been a significant issue in direct imaging of exoplanets, orbital motion has been 
considered in several closely related contexts. Here we briefly review a select portion of the literature. A very 
similar problem has been addressed in the context of searching for objects in our solar syste m, such as Kuiper 
Belt objects (KBOs), which can have proper motions on the order of 1" to 6" per hour (jChiang fc Brown 
19991 ). Blinking images to look for moving objects by eye is a well established technique. A more com- 
putationally intensive form of blinking images proceeds by shifting-and-adding a series of short exposures 
along trial paths, usually assumed to be linear. This "digital tracking" makes it p ossible to detect KBOs 
too faint to appear in a single exposure. This has been done both froni the ground dChiang fc: Brown 1999t 
Yamamoto et al.ll2008 ) and from space with HST ( Bernstein et al.ll2004 ). More recentlv lParker &: Kavelaard 
(|2010[ ) have taken into account nonlinear motion and optimized selection of the search space, especially 
important given the large data sets that facilities such as the Large Synoptic Survey Telescope will produce. 



O rbital motion is an important consideration when planning coronagraphic surveys of the HZs of nearby 



Brownl ( 20051 ) treats the problem of completeness extensively. Large parts of the HZ will be within 



stars. 

the inner working angle of the Terrestri al Plane t Find er-Coronagraph (TPF-C) and so undetectable during 
a single observation. Also discussed in iBrownl (J2005l) is photometric completeness - that is how long the 
TPF-C must integrate on a given star to detect an earth-like planet in the HZ. Other work on this topic 
includes iBrown &: Soummeii ( 20101 ) and iBrownl ( 20041 ) . These analyses consider orbital motion only between 
observations, not during a single observation as we do here. In general, the scenarios considered for these 
studies involved space-based high-performance coronagraphs on medium to large telescopes. In such cases 
exposure times were short enough and continuous so that orbital motion should be negligible during a single 
observation. 



The work most similar to our analysis here is the detection of Sirius B at lO/xm by Skemer &: Close! 

(|2011[ ). in fact, it was part of our motivation for the present study. ISkemer fc Close! (|201l!) used the well 



known orbit of the white dwarf companion to Sirius to de-orbit 4 years worth of images. Before accounting 
for orbital motion, Sirius B appeared as only a low SNR streak, but after shifting based on its orbit it 
appears as a higher SNR point source from which photometry can be extracted. Similar to this method, 
we will analyze the prospects for de-orbiting sequences of images, only we consider the case with no prior 
information at all, and with orbital elements with significant uncertainties. 



3. Quantifying The Problem 

In this section we will quantify the effects of orbital motion on an attempt to detect an exoplanet. Our 
first step will be to determine how fast planets move when projected on the focal plane of a telescope. Then 
we'll illustrate the impact this motion will have on the SNR and the statistical sensitivity of an observation. 



3.1. Basic Equations 

We begin by considering a focal plane detector working at a wavelength A in ^m. The FWHM of the 
PSF for a telescope of diameter D in m, neglecting the central obscuration, is 

FWHM = 0.2063— arcsec. (2) 

If we are observing a planet in a face-on circular (FOC) orbit with a semi-major axis of a in AU at distance 
d in pc, its angular separation will be a/d arcsec. At the focal plane the projected separation will then be 

p = 4.847^ in FWHM. (3) 

Act 

We note that it will occasionally be convenient to specify p in AU instead of FWHM. When it is not clear 
from the context we will use the notation pau to denote this. 



The orbital period is P = 365.25-\/a-^/M* days around a star of mass M* in Af©. In one period, the 
planet will move a distance equal to the circumference of its orbit, 27rp, so the speed of the motion in a FOC 
orbit will bqj 



-- - --- 1 ^) (T) (^) m) (^1 '" ^™" ^--' <^' 



In the general case, the equations of motion in the focal plane are 

esin(/) (cos(ri) cos(aj + /) — sin(i7) sin(w + /) cos(i)) 



1 

VFOC 



1- 

(1 + e cos(/)) (cos(rj) sin(w + /) + sin(f7) cos(a; + /) cos(i)) 



1 
vfoc 



e sin (/) (sin (f2) cos(a; + /) + cos(r2) sin(w + /) cos(i)) (5) 



l-e- 

(1 + e cos(/)) (sin(J7) sin(a; + /) — cos(17) cos(a; -|- /) cos(i)) 



Worn = vi^ijr^ 

where 51 is the longitude of the ascending node, lo is the argument of pericenter, i is the inclination, and 
the true anomaly / dep ends on a, e, and the time of pericenter passage r through Kepler's equation 



(JMurrav fc Correiall2010f ) 



In Figure[T]we show the variation in projected orbital speed for both circular orbits at several inclinations, 
and face-on eccentric orbits (i = 0), for a planet orbiting a IMq star at 1 AU. In the plots we normalized 
speed to 1, and provide vpoc for several interesting cases. These various scenarios produce projected orbital 
speeds of appreciable fractions of a FWHM per day. We will later show that, especially for ground based 
imaging, this causes a significant degradation in our sensitivity. 

Our main focus here is on planets in the HZ. Our simple definition of the HZ results in anz oc ^/IjI. 
Now, on the main sequence mass and luminosity approximately follow scaling laws of the form L* ex M^", 



*This result is equivalent to defining the gravitational constant in the focal plane as G = (0.0834Z)/(A(i))^ and using the 
equation for speed in a circular orbit Vdrc = \ -■ 
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Fig. 1. — Magnitude of projected orbital speed, normalized to 1 FWHM day ^, for 1 AU orbits around a IMq 
star. In (a) we show the orbital speeds for circular orbits at various inclinations, and in (b) we sho w the speeds 



for fa ce-on orbits at various eccen tricities. We give scali ng factors in (a) for M agAO/ VisAO (IClose et a.1. 



20121 ). GPI (iMacintosh et al.ll2012l). SP HERE/ZIMPOL (|Roelfsenia et al.ll2010l ). GMT (jJohns et al. 



201 



I. 



and E-ELT/EPICS (JKasper et al.ll2010l ). These scalings can be applied to the y-axis of either plot for various 
scenarios. These cases can also be scaled for different semi-major axes, telescopes, wavelengths, star masses 
and distances, by vpoc oc j^\/-^- See the text for the general equations of motion for arbitrarily oriented 
eccentric orbits. 



where 6 > 2 except for very massive stars. So according to Equation Q we expect vpoc in the HZ to 
increase as Af* decreases, i.e. M stars will have faster HZ planets than G stars. For example, a planet in the 
HZ of a Cen B (M* = O.9M0, L, — 0.5Lq) will be moving r oughly 20% faster t han a planet in the HZ of a 
Cen A (A/, = I.IA/0, i* = l.Si©) (stellar parameters from lBruntt et al.l ( 20101 )). 



To provide a more c oncrete example we return to the 20 hour observation of Gl 581d by the E- 



ELT/EPICS proposed bv iKasper et al.l (|2010l ). Using a wavelength of 0.75/iTO with Equation Q we find 
vfoc = 0.82 FWHM per day, or a total of 0.68 FWHM for a continuous 20 hour observation. Since this 
is a ground based observation the actual amount of motion to consider is ~ 1.15 FWHM over the ~ 1.4 
days minimum it wo uld take to integrate for 20 hours. Were this a face-on orbit, an eccentricity of 0.25 



(JForveille et al.ll201l[ ) would increase the maximum orbital speed to as much as 1.05 FWHM per day, or 1.47 



FWHM minimum for a 20 hour ground based observation. 



3.2. Impact on Signal-to-Noise Ratio 

So what does the orbital motion calculated above do to our observations? To find out we consider a 
simple model of aperture photometry. Let us assume that we are conducting aperture photometry with a 
fixed radius rap, that the PSF is Gaussian, and that we are limited by Poisson noise from a photon fiux N 
per unit area. With these assumptions, the optimum rap is 0.7 FWHM, but taking into account centroiding 
uncertainty rap ~ 1 FWHM is typical. We will approximate orbital motion at speed Vom by substituting 
a; — > a; — Vomt — xq. Orbits are of course not linear, but this will be approximately valid over short periods 
of time. The parameter xq allows us to optimize the placement of the aperture to obtain the maximum 
signal, i.e. centering the aperture in the planet's smeared out fiux. Note that with the exception of this 
centering parameter, this model appears quite naive in that we are not adapting the aperture radius and are 
pretending that we won't notice a smeared out streak in our images. 

Now the SNR in the fixed-size aperture after time At will be 

r^ap /'27r pAt 

SNRf,, = ■^" ^" ^° , (6) 

^NnrlpAUnt 

where Iq is the peak value of the PSF. In the case of no orbital motion Vom — and aperture rap — 1 FWHM, 
so we have 



SNRa='-^^. (7) 



As a simple alternative to a fixed size aperture, we also consider allowing our photometric aperture to expand 
along with the motion of the planet. This aperture will collect the same signal as in SNRo, but the noise 
increases with the area as 2rapVomAt, so we have 

SNR^^, = ^ »-^^°^ . (8) 

" ^N (1 + {2/7r)va^At) 

A convenient scaling is to multiply top and bottom by y/Vom and work in normalized SNR units of 
lo/y/Nvom- This puts time in terms of FWHM of motion, e = VomAt, and allows comparisons without 
specifying Vom- 



1.4 



1.2 



1.0 



5 0.8 



— 0.6 



0.4 



0.2 -/ 



0.0 



1.0 



0.9 



0.8 



of 

I 0.7 

01 



I 0.6 



0.5 ~- 



0.4 



0.3 



No Orbital Motion 
Orbital Motion, Fixed Ap. 
Orbital Motion, Expanding Ap. 



1.1' 



tx ^-1/z 



2 3 

Time (FWHM of motion) 



Orbital Motion, Fixed Ap. 
Orbital Motion, Expanding Ap . 



2 3 

Orbital Motion (FWHM) 



Fig. 2. — Top panel: SNR of a Gaussian PSF with and without orbital motion, in normalized units with 
time given as FWHM of motion. With no orbital motion SNRo ex ^/t. Equation ([6]) was used to calculate 
the SNR with orbital motion. After ^ 2 FWHM of movement, a maximum is reached and the observation 
can only be degraded by integrating further. Note that the fixed- aperture orbital motion case eventually 
goes down as SNR ex. 1/Vi- For comparison we also show the results with an aperture expanding with 
the moving planet, which eventually reaches a limit of 0.75. In the bottom panel we show the fractional 
reduction in SNR due to orbital motion for the fixed radius photometric aperture. 



In Figure [2] we plot the normalized SNR vs. time (measured in terms of FWHM of motion) with and 
without orbital motion and for both the fixed and expanding aperture cases. For the fixed aperture, after 
~ 2 FWHM of orbital motion a maximum of 0.69 is reached, and from there noise is added faster than 
signal. This means that further integration only degrades the observation. 

The expanding aperture SNRexp exceeds the maximum of SNRji^ after about 8 FWHM of motion, 
and 

lim 0.6 , v^ ^ ^ 0.6J- w 0.75. (9) 

^^°° Vl + (2/7r)a; V 2 

So if we integrate 4 times longer, adjusting the aperture size would allow us to gather a little more 
SNR, but only to a point. Given this large increase in telescope time for a relatively small improvement in 
SNR (only ^ 9% even if we integrate forever), and its better performance for smaller amounts of motion, 
the fixed-radius aperture will be our baseline for further analysis ~ keeping in mind that in some cases it 
may not be the true optimum. 

The peak in SNRfix (equation [5]) sets the maximum nominal integration time before orbital motion 
will prevent us from achieving the science goal. That is At,nax — {SNRmax/O-G)^- If the observation of a 
stationary planet would require an integration time longer than Atmax, then we can't achieve the desired 
SNR on an orbiting planet. This also sets the maximum orbital motion Smax = VomAtmax- From Figure [5] 
we find that emax — 1-3 FWHM. If more than 1.3 FWHM of motion occurs during an observation, we will 
not achieve the required SNR. 

We also show the fractional reduction in SNR in Figure [21 Almost no degradation occurs until after 
-- 0.2 FWHM of motion has occurred. SNR is reduced by -- 1% after 0.5 FWHM of motion, -- 5% after 
1.0 FWHM, and by -- 19% after 2.0 FWHM of motion. We must now decide how much SNR loss we can 
accept in our observation. 

The above analysis assumes a continuous integration. On a ground-based telescope one must consider 
that the maximum continuous integration time is < 12 hours, and in practice will likely be much shorter 
when performing high contrast AO corrected imaging. For instance, an exposure of 20 hours might have to 
be broken up over 4 or 5 or more nights, when considering the vagaries of seeing (required AO performance), 
airmass (either through transmission or rp requirements), rotation rate (for ADI), and weather. We can 
adapt the calculations for a ground based integration as follows 



'0 JO 
SNRnnd ~ 



J=*^ fti+At 



y^ / J g(-41n2((rcose-u„„t-a:o)^+r=sin2e))^^ 



j=i -h 



dOrdr. 



NnrlpAt, 



(10) 



In this expression we have broken the observation up into M integration sets which start at times tj and have 

lengths Atj . The total integration time is Atmt = /^ Atj and the total elapsed time of the observation is 

i=i 

Attot = tM + AtM — ti. 

We plot the results for a few ground-based scenarios in Figure |31 As one can see, observations of 
planets with orbital motion will be significantly degraded from the ground. This problem, which has been 
negligible in the high contrast planet searches to date, only becomes worse as we consider larger telescopes 
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Fig. 3. — Here we show the impact of orbital motion when combined with finite nightly integration times. 
The SNR of a Gaussian PSF with and without orbital motion is plotted in arbitrary units vs Atmt- The 
orbital speed Vo7n is given in FWHM day~^. 
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and improvements in AO technology which allow searches at shorter wavelengths. We next analyze how 
this reduction in SNR will affect our ability to detect exoplanets by increasing the rate at which spurious 
detections occur. 



3.3. Impact on Statistical Sensitivity 

Now we turn to the problem of detecting a planet of a given brightness. A planet is considered detected 
if its flux is above some threshold SNRt, which is chosen for statistical significance. The goal in choosing 
this threshold is to detect faint planets while minimizing the number of false alarms. For the purposes of 
this analysis we assume Gaussian statistics, in which case the false alarm probability (Pfa) per trial is 

Typically, planet hunters use a threshold of SNR — 5, which gives Pfa = 2.9 x 10"^. The number of false 
alarms per star, the false alarm rate (FAR), is then 

FAR = PFAXNt„als- (12) 



where Ntriais is the number of statistical trials per star. Following iMarois et al.l ()2008a|) . for a stationary 
planet Ntriais is just the number of photometric apertures in the image. A typical Nyquist sampled detector 
of size 1024x1024 pixels has Ntriais ~ 8 x 10''. Thus, an SNR = 5 threshold will result in FAR ~ 0.02 - 
about 1 false alarm for every 50 observatio ns. In the speckle lim ited case with non-Gaussian statistics, FAR 



will be worse than this for the same SNR (JMarois et al.ll2008a[) . In any case, the FAR is the statistic which 



determines the efficiency of a search for exoplanets with direct imaging. A high FAR will cause us to waste 
telescope time following up spurious detections, while raising the SNR threshold to counter this limits the 
number of real planets we will detect. 

The reduction of SNR caused by orbital motion confronts us with three options. Option I is to maintain 
the detection threshold constant and accept the loss of sensitivity. Option II is to lower the detection 
threshold to maintain sensitivity, accepting the increase in FAR. Option III is to correct for orbital motion, 
which as we will show also causes an increase in FAR. 



3.3.1. Option I: Do Nothing 

The default option is to do nothing, keeping our detection threshold set as if orbital motion is not 
significant. The drawback to this is that we will detect fewer planets. To quantify this we use the concept 
of completeness, that is the fraction of planets of a given brightness we detect. For Gaussian statistics and 
detection threshold SNRt = 5, the search completeness is given by 

n( \ 1 If f SNR{e)-5 \ 

Cie) = 1 - ^eric ( -= I . (13) 

where e = Vom^t is the amount of motion. In Figure S] (top) we show the impact of orbital motion on 
search completeness. Maintaining the detection threshold lowers completeness. How much depends on 
the completeness level, with brighter planets being less affected. For planets bright enough to yield 95% 
completeness with no motion, significant reduction in the number of detections begins after ~ 1 FWHM of 
motion. For 99.7% completeness the impact becomes significant after ^ 1.5 FWHM. 
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Fig. 4. — Top panel: completeness as a function of orbital motion if we maintain our detection threshold at 
5a. Planet brightness is expressed as the SNR at which we would be 50%(5cr), 68%(5.47(t), 95%(6.65(t), 
and 99.7% {7.75a) complete with no orbital motion. Bottom panel: the increase in false alarm probability 
(Pfa) if we lower the detection threshold to maintain 50% completeness for an orbiting planet that would 
have a brightness of 5a were it stationary. After ^ 1 FWHM of motion PpA increases exponentially until 
^ 4 FWHM where it becomes asymptotic to 0.5. 
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3.3.2. Option II: Lower Threshold 

Once orbital motion is recognized to be significant, a simple countermeasure would be to lower the 
detection SNR threshold in order to maintain completeness. The drawback to this option is that we have 
more false alarms, which must then be followed up using more telescope time. This results in a less efficient 
search. In Figure H] (bottom) we show PpA as a function of orbital motion, and denote the detection 
threshold we must use to maintain 50% completeness for a planet bright enough to give SNR — 5 were it 
stationary. Note that Pfa begins to increase exponentially after ^ 1 FWHM of motion. After ^ 4 FWHM 
PpA begins approaching 0.5 asymptotically. Once e sa 2 FWHM the number of false alarms per 1024x1024 
image approaches 1. 



3.3.3. Option III: De-orbit 

Option III is to correct for orbital motion, hoping to maintain sensitivity while limiting the increase in 
PpA- The essence of any such technique will be calculating the position of the planet during the observation, 
and de-orbiting in some way, say shift-and-add (SAA) on a sequence of images. The drawback of this 
approach is that it will produce more false alarms per observed star due to the increased number of trials, 
similar to lowering the detection threshold. If the orbit were precisely known, we could proceed with almost 
no impact on FAR. However, in the presence of uncertainties in orbital parameters or in a completely blind 
search we will have to consider many trial orbits. For now we can perform a "back-of-the-envelope" estimate 
of the number of possible orbits to understand how much FAR will increase. To do so, we begin by placing 
bounds on the problem. 

We can first establish where on the detector we must consider orbital motion. At any separation r from 
the star, the slowest un-bound orbit will have the escape velocity. Since we know that physical separation 
is greater than or equal to projected separation, r > p, and that maximum projected speed will occur for 
inclination i = 0, we know that 

Wesc = V2vFoc{a -^ p) (14) 

sets the upper limit on the projected focal plane speed of an object in a bound orbit. We can also set an 
upper limit on the amount of motion Cmax we can tolerate over the duration /Sttot of the observation based 
on the SNR degradation it would cause. So we only need consider orbital motion when 

V2VFOcip)^ttot > f^max- (15) 

From here we determine the upper limit on projected separation from the star for considering this problem: 

p™a. = 0.0136M, f^^^) AU. (16) 

By the same logic, for any point closer than pmax the maximum possible change in position is 

Apn.ax ~ V2VFOC (p) ^Uot m FWHM . (17) 

Then we must evaluate possible orbits ending anywhere in an area of 'K^/S.pmaxY' FWHM^ around an initial 
position. 

These two limits set the statistical sensitivity of an attempt to de-orbit an observation. The number of 
different orbits, Norb, will be determined by the area of the detector where orbital motion is non-negligible. 
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and the size of the region around each point that we consider. That is 

Pmai: 



Norbcc / Ap^^^pdp. (18) 







maWd-' 



e 



Xd 



Norb'^i — ] i^A A4,. (19) 



In general Ntriais oc No^b, so FAR ex Pfa x N^rb- Larger D, shorter A, closer d, and smaller acceptable 



\A 



orbital motion e will then all increase _F^il]J. Perhaps the most important feature of this result is that 
Norb (X Ai^Qj - increasing integration time rapidly increases the FAR of a blind search. Note that 
this is still less severe than the exponential increase in Pfa found for merely lowering the threshold. In the 
next section we will test these relationships after fully applying orbital mechanics, and see that they hold. 



4. Blind Search: Recovering SNR after Orbital Motion 

In this section we consider in detail a blind search, i.e. an observation of a star for which we have 
no prior knowledge of exoplanet orbits. We showed above that the problem is well constrained. Here we 
derive several ways to further limit the number of trial orbits we must consider. After that, we describe an 
algorithm for determining the orbital elements that must be considered and then discuss the results. Finally, 
we use this algorithm to de-orbit a sequence of simulated images and analyze the impact of correlations 
between trial orbits on FAR. 

To provide numerical illustrations throughout this section we consider the problem of a 20 hour obser- 
vation of a Cen A using the GMT at 10/im. This scenario is loosely based on performanc e predictions mad e 



for the proposed TIGER instrument, a mid-IR diffraction limited imager for the GMT (jHinz et al.ll2012l ). 
The details of these predictions are not important for our purposes, so we will only assert that this is a 
plausible case. There are other examples in the literature with similar integration times, such as the EPICS 
prediction we discussed earlier. 

We assume that this 20 hr observation is broken up into five At = 4 hr exposures, spread over 7 nights 
or /Sttot — 6.2 elapsed days from start to finish. The choice of Af is essentially arbitrary, but we have 
good reasons to expect it to be shorter than an entire night. An important consideration is the planned 
use of ADI, and the attendant need to obtain suffici ent field rotation in a short enough time to provide 



good PSF calibration while avoiding self-subtraction (jMarois et al.l 120061 ). The effect of airmass on seeing 



through ro oc cos(2:)^'^, where z is the zenith angle, and hence on AO system performance, could also cause 
us to observe as near transit as possible. Efficiency will be affected by chopping and nodding, necessary for 
background subtraction at 10/im. This will limit the net exposure time obtainable in one night.. 

Few ground-based astronomers would object to an assertion that we loose 2 nights out of 7 to weather. 
We could be observing in queue mode, such that these observations are only attempted when seeing is at 
least some minimal value, or precipitable water vapor is low. One can even imagine the opposite case at 
10/J,m, such that nights of the very best seeing are devoted to shorter wavelength programs. While this 
scenario may be somewhat contrived, we feel that it is both plausible and realistic. We now proceed to 



^Assuming background limited photometry with a diffraction hmited PSF, we expect Ai oc l/D* llHardvlll99a) . All else 
being equal, larger telescopes are better when considering this problem 
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describe a technique that would mitigate the effects of orbital motion for our GMT example and should be 
applicable to other long exposure cases. 



4.1. Limiting Trial Orbits 

Here we derive limits on the semi-major axis and eccentricity of trial orbits to consider. These limits 
are based only on the amount of orbital motion tolerable for the science case, and do not represent physical 
limits on possible orbits around the star. 

It is always true that r > p. This implies that, for any orbit, the separation of apocenter must obey 
Ta ^ P- This allows us to set a lower bound on a, Cmm, given a choice of e through 

Pau < aminC^ + e) (20) 

which gives 

armn = 0.2063—^^. (21) 

L)(l + e) 

The fastest speed in a bound planet's orbit will occur at pericenter, and using the maximum tolerable 
motion emax during our observation of total elapsed time Attot we can set an upper bound on a by noting 
that 

VFOcio-max)\ ^ttot ^ ^rnax (22) 

which leads to 

0.0834V ^^^*f^'- (23) 

XdJ 1 - e V ernax ) 

Using the GMT example: for e = 0.0, a^ax = 3.9 AU; and for e = 0.5, amax = 11-8 AU. Using Equation 
[T6]we have a projected separation limit of Pmax — 7.7 AU, so it is possible for these definitions to produce 
o-max < CLmin for Certain choices of e at a given p. This condition tells us that at such a value of e no orbits 
can move fast enough to warrant consideration. Thus we can set a lower limit on e at projected separation p 

e„^^n = \^/e+H-l-\ (24) 

where we have simplified by pulling out 

«-"«^(s)X^)'^ '^^' 

In practice, we might consider eccentricity ranges with Cmax less than 1, thus improving our sensitivity. 
Inputs to our choice of emax could include some prior distribution of eccentricities, or dynamical stability 
considerations in binary star systems and systems with known outer companions. 



4.2. Choosing Orbital Elements 

Now we describe an algorithm for sampling the possible trial orbits over a set of M sequential images. 
For now, we assume no prior knowledge of orbital parameters. We will employ a simple grid search through 
the parameter space bounded as described above. 
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1. Determine the region around the star to consider using Eq. (jl6[) . 

2. Identify regions of interest. In the best cases the orbital motion wih be smah enough that we wiU be 
able stack the images and search the result for regions with higher SNR (e.g. SNR > 4) and limit 
further analysis to those areas. In the worst cases orbital motion will be large enough that we will 
need to blindly apply this algorithm at each pixel within the bounding region identified in the previous 
step. In the present GMT-aCen example we are in the former case. 

3. For each region, choose a size, perhaps based on v^sc (as in Eq. [T7|) . 



4. Chose a starting point {xi,yi), with pi — \/ x\ + y\. If we are proceeding pixel by pixel, then (xi,2/i) 
describes the current pixel. 

5. Choose e € eminipi) ■ ■ ■ Smax using Equation (p4|) and assumptions about Cmax- 

6. Choose a € amm(pi,e) . . -amaxie) using Equations (|211) and ([23)1 . 

7. Choose time of pericenter r e ii — P{M^,,a) . . .ti where P is the orbital period and ti is the time of 
the first image. Now calculate the true anomaly /(ii; a,e,r, P) using Kepler's equation and physical 
separation using: 

a(l-e) 



1 + ecos(/) 



(26) 



8. if e 7^ 0: Choose w e . . . 27r 
if e = 0: set a; = 0. 

9. if sin(a; + /) > 0: 

(a) Given e, a, t, /, and cu, calculate 



±A 



sm(a; + /) 



Sm i Z — , ^, r\ . ■ 2/ \ — ^^ , ,9 -x (,^°J 



COS Si = 7-5 TT— (29) 



y cos(w + f) — X sin(a; + /) cos i 
r(cos^(a; + /) + sin^(aj + /) cos^ i) 

y sm{Ld + /) cos i + x cos(a; + /) 
r(cos^(aj + /) + sin^(a; + /) cos-^ i) 
where f2 should be determined in the correct quadrant. 

(b) We now have a complete set of elements, and so can SAA the sequence of images based on these 
orbits (one for each i). Doing so requires calculating the true anomaly fj at the time of each 
image, and then calculating the projected orbital position of the prospective companion in each 
image. 

10. if sin(w + /) = 0, we do not have a unique solution for inclination. This is the special case where the 
planet is passing through the plane of the sky. 
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or for w + / = TT calculate ft: 



(a) for w + / = calculate J7: 

(30) 

(31) 

(32) 

(33) 

determining il in the correct quadrant. 

(b) Choose i £ . . .n 

(c) We now have a complete set of elements, and so can SAA as in step I9bl above. 

(d) Repeat steps llObI to llOcI until all i chosen. 

11. Repeat the above steps until the parameters oj, r, a, and e are sufficiently sampled for each starting 
point. 



sinfi = 


y 
r 


cosfl = 


X 

r 


sinO = 


-y 

r 


cosfJ — 


~x 



4.3. De-orbiting: Unique Sequences of Whole-Pixel Shifts 

The algorithm just described will produce a large number of trial orbits, many of which will be very 
similar. The information content of our image is set by the resolution of the telescope, so we can take 
advantage of this similarity to greatly reduce the number of statistical trials. This is done by grouping similar 
orbits into sequences of whole-pixel shift sequences, where the pixels are at least as small as FWHM/2. As 
we will see, we typically will want to oversample, to say FWHM/3, to ensure adequate SNR recovery. 

We calculate the pixel-shift sequence for each orbit by determining which pixel the trial planet (or 
rather, the center of its PSF) lands on at each time step. Many orbits end up producing the same sequences 
of pixel-shifts, and we will keep only the unique ones for use in de-orbiting the observation. In Figure [S] we 
illustrate the outcome of the pixel-shift algorithm, showing two unique sequences and a few of the orbits 
that produced them. 

To test the above algorithm and the pixel-shift technique, we used our GMT a Cen A example and 
determined the trial orbits for various separations and Ais. We set Smax =0.5 based on our earlier analysis 
of SNR. The results are summarized in Figure |6l The problem is generally well constrained in that we only 
have a finite search space for any initial point. The data used to construct Figure [S] are provided in Table 
[TJ Comparing Norb to Nghifts, note the large reduction in the number of trials (~ 10^ to ^ 10^) due to 
combining similar orbits. 



4.4. Norb Scalings 

In Figure [7] we plot the area of the detector which contains the possible trial orbits at pi = 1.0 AU 
vs. the total elapsed time Attot- We conclude from this plot that the area around a given starting point 
is proportional to At^^j. Also in Figure [7] we plot area vs separation from the star, and conclude that area 
is proportional to l/pi- Taken together these results give confidence that the N^rb fx Atf^t scaling derived 
earlier holds when we fully apply orbital mechanics rather than the escape velocity approximation. 
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Fig. 5. — Two sequences of whole pixel-shifts, one in red and one in blue. We also show a few of the 
many orbits that produce these shift sequences. Once these shifts are determined, a set of 5 images can be 
de-orbited by shifting the images by the indicated sequence 5-4-3-2-1, that is the pixel containing the orbit 
in image 2 is is shifted and added to the pixel containing the orbit in image 1, and likewise for images 3, 4, 
and 5. Of course, the entire image is shifted, not just single pixels. 
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Fig. 6. — Example trial orbits for the GMT, working at lO/im, observing aCen A. Plotted are the end points 
of orbits calculated using the algorithm given in Section F4.2l for the given initial projected planet separations 
pi and elapsed observation times A.ttot- The red points show the effect of changing initial separation for 
a constant elapsed time. At 1 AU initial separation the colors correspond to different elapsed times as 
indicated in the legend. We further analyze these relationships in Table [T] and Figure [T] The results of 
the algorithm appear more complicated than the simple escape- velocity circle analysis in Section [3.3.3l The 
end-point clouds are not circularly symmetric about the starting point, and have some azimuthal structure. 
For instance there is a triangle extending azimuthally corresponding to face-on high-e orbits, and there are 
gaps along the radius from the star corresponding to i very near 90°. These structures are consequences of 
the chosen grid resolution. 
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Pi (AU) 


Attot (days) 


No. Obs. 


Norb 


Nshifts 


0.5 


6.0 


5 


2.7 X 10** 


285 


1.0 


2.0 


5 


4.1 X 10^ 


14 


1.0 


4.0 


5 


4.1 X 10^ 


76 


1.0 


6.0 


5 


4.1 X 10^ 


134 


1.0 


8.0 


5 


4.1 X 10^ 


253 


1.5 


6.0 


5 


5.2 X 10^ 


90 


1.0 


2.0 


3 


4.1 X 10^ 


10 


1.0 


4.0 


5 


4.1 X 10** 


78 


1.0 


6.0 


7 


4.1 X 10** 


292 


1.0 


8.0 


9 


4.1 X 10^ 


815 



Table 1: Results of applying the algorithm detailed in Section |4~2] for various separations and elapsed obser- 
vation times. See also Figure[7l Note the dramatic reduction in the number of trials {Norb vs. Nshifts) after 
combining similar orbits into whole-pixel shift sequences. 
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Things are a bit more complicated when we consider the scahng of the number of while-pixel shift 
sequences. We conducted two sets of trials at pi — 1.0 AU. In the first, the number of observations and 
their relative spacing was held constant regardless of Aitot. In the second set, the number of observations 
scaled with Atfot. As shown in Figure [71 when the number of observations is constant, the number of shifts 
scales as Ai^ ^^ , but when the number of observations grows with Aijot the number of shifts scales as roughly 
AijQj. Figure [Tjl shows that the number of shifts scales as 1/pi. Taken together, we see that for a constant 
number of observations the pixel-shift technique will follow the A'orb oc Atj^j scaling. However, if the number 
of observations also scales with tdtot , then our results imply that Norb oc At^^f . The value of the exponent 
likely depends on the details of the observation sequence, but this has important implications for observation 
planning. 



4.5. Recovering SNR 

We next consider whether de-orbiting by whole-pixels adequately recovers SNR. To test this we "or- 
bited" a Gaussian PSF on face-on orbits with various eccentricities, starting from pericenter. We then 
calculated shifts for detector samplings of 2, 3, and 4 pixels/FWHM, and then de-orbited by these shifts. 
The results are summarized in Table [2] On a critically sampled detector we only recover a 5a planet to 
~ 4.9a, a 2% loss of SNR. At 3 pixels/FWHM we do much better, recovering SNR to 4.97 for low eccen- 
tricities, and 4.95 for higher eccentricities. Performance for 4 pixels/FHWM sampling is similar. A 2% loss 
of SNR nearly doubles Pfa, so it appears that we should oversample to at least 3 pixels/FWHM, either 
optically or by re-sampling images during data reduction. In our analysis we have assumed that the limiting 
noise source is background photons (PSF halo or sky), so we ignore the increased readout noise expected 
from oversampling. 



4.6. Correlations And The True Impact On Ppj^ 

As we have noted several times, the main impact of orbital motion is to reduce SNR, which in turn 
reduces our statistical sensitivity. If we attempt to de-orbit an observation in order to recover SNR, we 
do so at the cost of a large increase in the number of trials. Worst case, this results in a proportional 
increase in FAR since nominally FAR = Pfa x Norb- However, we expect significant correlation between 
trials of neighboring orbits and whole-pixel shifts. To investigate this, we performed a series of monte carlo 
experiments. A sequence of images with Gaussian noise was generated, and first stacked without shifting, 
hereafter called the naive-add. The same sequence was then shifted by each possible whole-pixel shift. 

Table 2. SNR recovered after de-orbiting with whole-pixel shifts for various 

samplings. 



Sampling 






SNR Recovered 






(pix/FWHM) 


e=0.0 


e=0.1 


e=0.2 e=0.3 e=0.5 


e=0.7 


e=0.9 


2 


4.89 


4.89 


4.89 4.88 4.86 


4.86 


4.86 


3 


4.97 


4.96 


4.95 4.94 4.94 


4.95 


4.95 


4 


4.97 


4.97 


4.97 4.97 4.94 


4.92 


4.92 
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Fig. 7. — Scaling of the number of orbits and the number of resuhing whole-pixel shifts with observation 
elapsed time and with distance from the star. These results demonstrate that the number of trial orbits 
Norb oc Ai^oj scaling that we derived using the escape velocity holds when we rigorously apply orbital 
mechanics. Note though that the situation is more complicated with the number of shifts - if the number 
of observations increases with elapsed time then the number of shifts grows faster than Afj^^j, implying that 
Norb will increase faster than Atf^t- These scalings lead to one of our main, if seemingly obvious, conclusions: 
one must limit the elapsed time of an observation as much as possible when orbital motion is significant. 
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assuming a lAU initial separation around a Cen A. This experiment was conducted for observations with 
total elapsed times Attot of 4.2, 6.2, 8.2, and 10.2 days, with samplings of 2, 3, and 4 pixels/FWHM. 

We performed several tests on each sequence. The first was a simple threshold test on the naive-add, with 
the threshold set for the worst case orbital motion given by Equation 1101 with Vom — Vesc- We performed 
simple aperture photometry, with a rap = 1 FHWM. As expected the resultant Pfai is as predicted by 
Equation [TTJ The next test was to apply a 5(t threshold after de-orbiting by whole-pixel shifts and adding. 
If all shifts were completely uncorrelated, then we would expect FAR — (2.9 x 10~^) x Ng^ifts, but as we 
predicted, shifts are correlated and PfA2 is lower than this. 

The final test performed was to apply both thresholds in sequence, such that a detection is made only if 
the naive-add results in SNR greater than the threshold for worst case orbital motion, and the de-orbited 
SAA results in SNR > 5. This PfA3 is lower than either Pfai or PfA2, but still higher than if no orbital 
motion occurred. 

The results of each trial are present in Table [3l Applying both threshold tests results in significant 
improvement over the naive-add in terms of FAR. Another interesting result is that sampling has only a 
minor impact on PfA3- This makes some sense as we expect the correlation of neighboring shifts to be set 
by the FWHM, not the sampling. So even though the accuracy of SNR recovery is improved, and quite a 
few more shifts are required, these shifts remain correlated across the same spatial scale resulting in little 
change in the overall FAR. 



4.7. Impact on Completeness of the Double Test 

There is still an impact on completeness, however, because we are now conducting two trials instead 
of one. This lowers the true positive probability (Ptp). Consider a 5a planet on the worst case fastest 
possible orbit, for the 10.2 day elapsed time case. The threshold for the naive add is 2.625. We have a 50% 
probability of detecting this planet after the naive add. If it is detected on the first test, there is then some 
probability Ptp < 1 of detecting at SNR > 5 after de-orbiting. Worst case, this will be 50%, resulting in a 
net Ptp of 25%. In reality, it will be better than this as the two trials will be strongly correlated. 

Even if this worst case of 25% were realized this is still significant improvement over Option I. A 2.6cr 
signal would only be detected 10% of the time with a 5a threshold. Given the reduction in PpA from 
4.3 X 10^'^ to 2.2 X 10~^, likewise an improvement over Option II at 2.6cr, it is clear that de-orbiting by 
whole-pixel shifts does improve our ability to detect an orbiting planet. The situation will be even better for 
slower planets, and most of the area searched will not be subject to the worst case orbital speed. We leave 
a complete analysis of the impact on search completeness for future work. One can also imagine adjusting 
the thresholds to optimize completeness at the expense of worse Pfa- 



4.8. Tractability of a Blind Search 

We end this section by concluding that a blind search when orbital motion is significant is tractable. 
Orbital motion will make such a search less sensitive, both in terms of number of false alarms and in terms of 
completeness, but Keplerian mechanics gives us enough tools to bound the problem. As we have shown de- 
orbiting a sequence of observations can recover SNR to its nominal value, and we can do so while controlling 
the impact on statistical sensitivity. For the Attot — 6.2 day observation, PfA3 was roughly a factor of 10 
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Table 3. False alarm probabilities after de-orbiting Gaussian noise images. 



Attot (days)i 


SNRt^ 


N.Mfts^ 


Pfai^ 


PfA2^ 


Pfa^' 




2 pixels/FWHM 




4.2 


4.635 


64 


1.74 X 


10-6 


7.65 X 


10-6 


8.06 X 10-^ 


6.2 


4.220 


122 


1.24 X 


10-5 


1.52 X 


10-5 


2.70 X 10-6 


8.2 


3.330 


231 


4.40 X 


10-4 


2.71 X 


10-5 


9.93 X 10-6 


10.2 


2.625 


364 


4.33 X 


10-3 


4.03 X 


10-5 


2.17 X 10-5 




3 pixels/FWHM 




4.2 


4.635 


108 


2.11 X 


10-6 


1.37 X 


10-5 


4.80 X 10-^ 


6.2 


4.220 


285 


1.21 X 


10-5 


3.39 X 


10-5 


2.04 X 10-6 


8.2 


3.330 


496 


4.31 X 


10-4 


5.64 X 


10-5 


9.96 X 10-6 


10.2 


2.625 


741 


4.34 X 


10-3 


8.15 X 


10-5 


2.67 X 10-5 




4 pixels/FWHM 




4.2 


4.635 


217 


1.78 X 


10-6 


2.64 X 


10-5 


4.44 X 10-^ 


6.2 


4.220 


487 


1.24 X 


10-5 


5.61 X 


10-5 


1.48 X 10-6 


8.2 


3.330 


844 


4.35 X 


10-4 


9.19 X 


10-5 


1.14 X 10-5 


10.2 


2.625 


1315 


4.32 X 


10-3 


1.41 X 


10-4 


3.15 X 10-5 



^Elapsed time of the observation. 

^SNR threshold from Equation [TUl using Vorb = V^vpoc- 

^Number of unique whole-pixel shifts required to de-orbit. 

4False alarm probability for the naive-add, from MC experiment results. Ex- 
pected values given by Equation [TT] 

5False alarm probability after de-orbiting with whole-pixel shifts. 

6False alarm probability after testing both the naive-add and de-orbiting. 
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higher than if no orbital motion occurred. This increase only occurs over a bounded region around the star, 
so the net effect on FAR will be contained. Using this factor of 10 as the mean value over the 7.7 AU = 69.1 
FWHM radius region around a Cen A where orbital motion is significant, the FAR in this area will have 
gone from ^ 1/1000 to ^ 1/100 in our GMT/10/im example. The key, though, appears to be to limit the 
elapsed time of the observation as the number of trials increases — decreasing sensitivity - proportionally 
to at least Atf^^ in a blind search. 

The main caveat at this point in our analysis is that we have drawn the conclusion of tractability using 
Gaussian statistics. It is well known that speckle noise, which will often be the limiting noise s ource for high 



contra st imaging in the HZ, is not Gaussian and results in much higher PpA for a given SNR (jMarois et al. 
2008al ). Future work on this problem will need to take this into account. 



Next we consider a more strongly bounded scenario, where we have significant prior information about 
the orbit of the planet from radial velocity surveys. 



5. Cued Search: Using RV Priors 

The situation is greatly improved if we have prior information, such as orbit parameters from RV or 
astrometry. Here we consider th e case of Gliese 581d , and the previously discussed future observation of this 
planet by EPICS at the E-ELT (JKasper et al.ll2010). There is some controversy surrounding the solut ion to 
the RV signal, and whether planet d even exists (IForveille et al. 2011 : Vogt et al. [20121 : lBaluevll2012[ ). We 



show results for both the floati ng eccentric i ty Ke plerian fits of IForveille et al.l (|201lh [hereafter Fll], and the 
all circular interacting model of lVogt et al.l (|2012l ) [hereafter VI 2]. Doing so allows us to illustrate the impact 
of eccentricity on the analysis, and prevents us having to take a stand in a currently raging debate. The 
parameters used herein are listed in Table 21 

Instead of a grid search, we use a monte carlo (MC) method. The RV technique provides the parameters 
a, e, Lo and io or their equivalents. We can take the results of fitting orbits to the RV signal, and the 
associated uncertainties, as prior distributions which we sample to form trial orbits. We will assume that all 
uncertainties are uncorrelated and are from Gaussian distributions. 

We assu me that the 20 hr in tegration is broken up over 6.2 nights based on the same logic discussed 
in Section m iKasper et al.l ( 20101 ) actually assumed 20 x 1 hr observations based on the amount of rotation 
needed, but did not consider the effects of orbital motion over 20 days of a 67 day period (M. Kasper, 
personal communication (2012)). 



Model 


a (AU) 


e 


w (deg) 


cTto (days) 


Forveille et al. (2011) 
Vogt et al. f2012) 


0.218 ±0.005 
0.218 ±0.005 


0.25 ±0.09 
0.0 ±0.0 


356.0 ±19.0 
0.0 ±0.0 


±3.4 

±7.45 



Table 4: Orbital parameters for Gl 581d used in this analysis. We derived the values reported here from 
other parameters where necessary. Only the uncertainty in t^ impacts our analysis. In both models the 
orbital period is 66.6 days. 
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5.1. Constraints 

In order to minimize the number of trial orbits to consider, we can apply various constraints taking 
advantage of the information we have from the RV detection. 

In the case of a multi-plane t system dynamical analysis can place constraints on the inclination based 
on system stability. For Gl 581, iMavor et al.l ( 20091) found the system was stable for i > 30. We can also 



make use of the geometric prior for inclination, where we expect Pi — sin(i) in a population of randomly 
oriented systems. 

Since this is a reflected light observation, the orbital phase and its impact on the brightness of the planet 
must be considered. The planet's reflected flux is given by 

F,{a) = F.(^^^ M\Ma) (34) 

where F^ is the stellar flux, Rp is the planet's radius, r its separation, Ag{X) is the wavelength dependent 
geometric albedo, and $ is the phase function at phase angle a. The phase angle is given by 

cos(a) = sin(/ + cj) sin(z). (35) 

In general, determining the quantity Ag{\)^{a) requires atmospheric modeling ( Cahov et al.l 120101) . For 
now, we assume that $ follows the Lambert phase function 

$(q;) — — [sin(a) + (tt — a) cos(q;)] (36) 

TT 



We assume that the prediction of iKasper et al.l (|2010l) was made for the planet at quadrature, a — tt/2, 



where $ = 0.318. We then require that the mean value of $ during the observation be greater than this 
value - that is the planet is as bright or brighter than it is at quadrature. 



5.2. Initial Detection 

An important consideration in an RV-cued observation will be when to begin. As a first approximation, 
we assume that maximizing planet-star separation will maximize our sensitivity. This may not be true when 
working in reflected light due to the phase and separation dependent brightness of the planet in this regime. 
Proceeding with the approximation for now, we expect to plan this observation to be as close to apocenter 
as possible. In this case we will begin integrating 3.1 days before to + P/2. 

To understand the area where we will be searching for Gl 581d, we first conducted an MC experiment 
to calculate the possible positions of the planet at t = to + -P/2 — 3.1 days. To do so, we drew random values 
of a, e, w, and to from Gaussian distributions with the parameters of Table |4l We drew a random value of i 
from the sin(i) distribution, and rejected any value of i < 30 based on the dynamical prior. Finally fl was 
drawn from a uniform distribution in . . . 27r. This process was repeated 10^ times, and the frequency at 
which starting points occur in the area around the star was recorded. The results are shown in Figure H] for 
the V12 circular model and for the Fll eccentric model. The Figure shows the area which must be searched 
to obtain various completeness. For instance, if we desire 95% completeness in the V12 model, we must 
consider an area of 71 apertures. Since this SNR = 5 detection is broken up into 5 distinct integrations, our 
first attempt will have SNR — 2.24, giving a FAR = 0.89 for the first 4 hr integration. In other words, we 
should expect a false alarm in addition to a real detection. 
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Fig. 8. — Pos s ible st arting points for Gl 581d, observed near apocenter. Top : using the param eters of 



Forveille et al.l (J2011r )'s eccentric model. Bottom: assuming the parameters of IVogt et al.l (|2012n 's circu- 
lar interacting model. The color shading is in units of probability per aperture (each aperture has area 
ttFWHM^). The legend indicates the color which encloses the given completeness intervals, and the enclosed 
area in apertures, which can be directly related to the false alarm rate as discussed in the text. 
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5.3. Calculating Orbits and Shifts 

Now we assume that we have an initial detection at SNR ^2.24 within the highest probability regiont|f|. 
In order to follow-up this detection over subsequent nights, we must determine the possible locations of the 
planet, constrained by the RV-derived orbital elements. 

We proceed by choosing a, e, w, and to from Gaussian distributions as above. Now as long as r > p 
we will have a unique solution for i and il given the randomly chosen parameters (see the blind search 
algorithm above). We take into account dynamical stability by rejecting any orbit which has i < 30. The 
orbit determined in this fashion was then projected 6.2 days into the future and the frequency of these final 
points was recorded. We show the result for the VI 1 model in Figure |9l top panel. Using the RV determined 
parameters and their uncertainties allows us to determine the probability density of orbit endpoints, and 
determine how much of the search space we must consider for a given completeness. The whole-pixel shifts 
were also calculated using a sampling of FWHM/3, and are shown in the legend. We also applied the blind 
search algorithm to this observation from the same starting point, and show the results for comparison in 
the bottom panel of Figure O As expected the RV priors significantly reduce the search space - we have 942 
trial shift-sequences to consider instead of 12000. 

Another important consideration here is that our initial 2.24(T detection will have a large position 
uncertainty, which we estimate by ap„ — FWHM/ SNR. We added a random draw for the starting position, 
and repeated the MC experiment for Fll and also conducted a run for the V12 parameters. The results are 
shown in Figure 1101 The number of shift sequences is much higher due to the uncertainty in the starting 
position caused by our low SNR initial detection, but we expect correlations to come to the rescue as in our 
a Cen example. To compare to Figure |9] keep in mind that the blind search would have to be applied to all 
5500 pixels in the search space indicated by Figure ID 

As in the GMT/a Cen example, we leave for future work a complete analysis of sensitivity and com- 
pleteness. The large number of trial shifts calculated when we include uncertainty in the starting position 
motivates us to suggest that we will ultimately turn this analysis over to a much more robust optimization 
strategy, such as a Markov Chain Monte Carlo (MCMC) routine. Once an area of the image was identi- 
fied with a high post-shift SNR, a MCMC analysis could determine the very best orbit and assign robust 
measures of significance to the result. 

We also note that these results likely overestimate the number of trial orbits since we have assumed 
uncorrelated errors. In reality the RV best fit parameters are likely strongly correlated, which should act to 
reduce the number of orbits to consider. 



6. Conclusions 

In the coming campaigns to directly image planets in the HZs of nearby stars, orbital motion will be 
large enough to degrade our sensitivity. This effect has been ignorable in direct imaging campaigns to date, 
which have typically looked for wide separation planets. We have analyzed this issue in some detail, and 
shown that applying basic Keplerian orbital mechanics allows us to bound the problem sufficiently that we 
believe direct imaging in the HZ to be a tractable problem. Our main conclusions are: 



-fPor the purposes of this analysis, we calculated initial separation pi using the mean parameters for each model and an 
inclination i = 60 
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Fi g. 9. — Trial orbits for Gl 581d, observed near maximum elongation. In the top panel we use the parameters 
of iForveille et al.l (|201ll )'s eccentric model. The bottom panel shows the results for a blind search from the 
same starting point. The red cross shows the starting point, and the star is located at the origin. The top 
panel color shading is in units of probability per aperture (each aperture has area ttFWHM^). The legend 
indicates the color which encloses the given completeness intervals, the enclosed area in apertures, and the 
number of unique whole-pixel shift sequences which must be tried in order to de-orbit the observation. 
The number of shift sequences is directly related to the false alarm rate, and hence the sensitivity. For 
comparison, the blind search algorithm produced ~ 12000 shifts. RV cueing greatly improves our sensitivity 
in the presence of orbital motion. 
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Fig. 10. — Trial or bi ts for Gl 581d, observed near maximum elongati on, assuming the parameters of 
(top) iForveille et al.l (|201l[ )'s Keplerian eccentric model and (bottom) IVogt et al.l (|2012l )'s circular in- 
teracting model. In this simulation we allowed the initial position to vary with standard deviation 
o'x,y — FWHM/SNR. The red cross shows the starting point, and the star is located at the origin. 
The color shading is in units of probability per aperture (each aperture has area ttFWHM^). The legend 
indicates the color which encloses the given completeness intervals, the enclosed area in apertures, and the 
number of unique whole-pixel shift sequences which must be tried in order to de-orbit the observation. The 
number of shift sequences is directly related to the false alarm rate, and hence the sensitivity. 
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(1) When projected onto the focal plane, a planet in a face-on circular orbit moves with speed given by 

(37) 



VFOC = 0.0834 (— j J^- FWHM day"\ 



In the HZ of nearby stars, especially when considering giant telescopes, speeds are high enough that planets 
will move significant fractions of a PSF FWHM during a single observation. This smears out the planet's 
flux resulting in a lower SNR. 

(2) In background limited photometry, an SNR maximum is reached after about ^ 2 FWHM of motion 
has occurred on the focal plane. From there, integrating longer offers no improvement with a flxed-size 
aperture. Adapting the aperture could mitigate this to some extent, but at the cost of significantly longer 
exposure times. 

(3) When SNR is reduced by orbital motion, we have three options. Option I is to do nothing, and 
accept the loss of completeness due to planets appearing fainter. Option II is to adjust our detection threshold 
at the cost of more false alarm detections. Option HI is to de-orbit an observation, recovering SNR to its 
nominal value, but also at the cost of more false alarms. 

(4) For exposure times of 10s of hours, we expect an observation to extend over several days under 
realistic assumptions about ground based observing. If we naively attempt to de-orbit such an observation, 
the false alarm rate per star will increase by at least FAR ex At^^j, where Attot is the total elapsed time of 
the observation. 

(5) De-orbiting a sequence of shorter exposures is possible, and tractable. Taking advantage of strong 
correlations between trial orbits, we will realize increases in the FAR on the order of a factor of 10 in the 
region around a star where orbital motion matters. Since this will be a small, bounded region, this increase 
in FAR appears to be acceptable. 

(6) Cueing from another detection method, such as RV, provides significant benefit. It allows us to 
initiate our search at the optimum time, and significantly reduces the size of the search space. Having prior 
distributions for some of the orbital elements will allow us to efficiently determine where and how to search 
to optimize completeness. 
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for reviewing this manuscript. JRM is grateful for the generous support of the Phoenix ARCS foundation. 
LMC and JRM acknowledge support from the NSF AAG. 
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